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Abstract. We calculate real-space static correlation functions related to basic 
entities of the one-dimensional Hubbard model, which emerge from the exact Bethe- 
ansatz solution. These entities involve complex rearrangements of the original 
electrons. Basic ingredients are operators related to unoccupied, singly occupied with 
spin up or spin down and doubly occupied sites. The spatial decay of their correlation 
functions is determined using an approximate mean-field-like approach based on the 
Zou-Anderson transformation and DMRG results for the half-filled case. The nature 
and spatial extent of the correlations between two sites on the Hubbard chain is studied 
using the eigenstates and eigenvalues of the two-site reduced density matrix. 



PACS numbers: 71.10.Fd, 71.10.Pm, 73.90.+f 
1. Introduction 

It was recently shown in [Ij that a consistent description of the Bethe-ansatz exact 
eigenstates of the Hubbard chain can be achieved in terms of some basic entities called c 
fermions, spinons and ry-spinons. (In the related preliminary studies of [2j such objects 
were named c pseudoparticles, spinons and holons, respectively.) These objects are 
related to the rotated-electron occupancy configurations. The rotated-electron creation 
and annihilation operators are related to those of the original electrons by a unitary 
transformation, V{U). This transformation is defined such that double occupancy of 
these rotated electrons is a good quantum number for all values of U. There are infinite 
choices for such transformations. Examples are those reported in [HJl^. However, the 
BA solution performs a specific electron - rotated-electron unitary transformation [1]. 

The operator formulation introduced in P accounts for all representations of the 
model global [S0{4) ® t/(l)]/^2 symmetry algebra. Here S0{4) = [SU{2) ® SU{2)]/Z2 
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refers to the spin SU{2) and ry-spin SU{2) symmetries and U{1) to the c hidden U{1) 
symmetry found in [5j . We denote the spin and the ry-spin of an energy eigenstate by Sg 
and Sj^^ respectively. We caU 2Sc the number of rotated-electron singly occupied sites, 
which is the eigenvalue of the generator of the c hidden U{1) symmetry. 

At very large U double occupancy is a good quantum number, since there is a 
very large energy separation between states differing by the number of doubly occupied 
sites (in this limit the unitary transformation is the identity) . As the Coulomb repulsion 
becomes finite, the unitary transformation is not known explicitly but it can be expressed 
in a perturbative expansion in powers of t/U ^ as shown for instance in [6j. To leading 
order mt/U the unitary operators associated with such transformations have a universal 
form. Such operators diflFer in their higher-order terms. Recently, the matrix elements 
of the unitary operator associated with the specific transformation performed by the 
Bethe-ansatz solution has been obtained for the whole U/t > range in the basis of the 
energy eigenstates [Ij. 

The various entities referred to above are introduced so that the number of c 
fermions Nc = 2Sc and that of c fermion holes A^^ = [Na — 2Sc\ are equal to the 
number of singly occupied sites of rotated electrons (with spin up or spin down) and 
the number of rotated-electron doubly and unnocupied sites, respectively. They store 
information on the charge part of these rotated electrons. The number of spinous 
Ms = 2Sc also equals that of singly occupied sites of rotated electrons. The spin-1/2 
spinous of component +1/2 have information on the spin of these rotated electrons 
of spin up, the spin-1/2 spinous of component —1/2 is related to the spin part of such 
rotated electrons with spin down. The number of ry-spinons = [Na — 2Sc] equals that 
of rotated-electron doubly and unnocupied sites. While the = [Na — 2Sc\ c fermion 
holes describe the c hidden U{1) symmetry degrees of freedom of such sites occupancies, 
the = [Na — 2Sc] ry-spinons refer to the 77-spin SU (2) symmetry degrees of freedom 
of the same site occupancies. Specifically, the 7y-spin-l/2 ry-spinons of component +1/2 
are related to the unoccupied occupied sites and the r7-spin-l/2 ry-spinons of component 
— 1/2 are related to the doubly sites of the rotated electrons. 

Furthermore in [K [2] it was proposed that both the energy eigenstates inside 
and outside the Bethe-ansatz solution subspace can be generated by occupancy 
configurations of the c fermions (associated with the usual real charge rapidities of 
the Bet he ansatz solution) and suitable spinon and ry-spinon occupancy configurations. 

Out of the Ms = 2Sc spinous, a number M^^ = 2Ss of spinous are unbound and 
determine the energy eigenstate spin value Ss- For the energy eigenstates inside the 
Bethe-ansatz subspace, all unbound spinous have spin projection +1/2. Flipping such 
unbound spinous generates the spin towers of spin projection = —Ss^—Ss + l^...^Ss — 
Ss- The corresponding energy eigenstates are outside the Bethe-ansatz subspace. The 
remaining M^^ = [2Sc — 2Ss] spinous are bound within composite entities with total spin 
zero (bound-state between a spin-1/2 spinon with component 1/2 and another spin-1/2 
spinon with component —1/2). There are also Bethe-anstaz excited states containing 
z/ = 1,2,3,... spin-neutral pairs of such bound spinous. Note that zero-magnetization 
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ground states have no unbound spinous, so that aU their spinous are bound within 
spin-neutral pairs. In this case one has that Mg = 2Nsi^ where A^^i is the number of 
spin-singlet two-spinon composite objects called in [Ij si fermions. 

Similarly, out of the = [A^^ — SS'c] ry-spinons, a number M^^ = 25^^ of ij- 
spinons are unbound and determine the energy eigenstate ry-spin value Sj^. For the 
energy eigenstates inside the Bethe-ansatz subspace, all unbound ry-spinons have ry-spin 
projection +1/2. Those correspond to rotated-electron unoccupied sites. Flipping such 
unbound ry-spinons generates the ry-spin towers of ry-spin projection = — /S^, — S'^ + 
1, Srq 1, Srq. Such ry-spin flipping processes involve creation of on-site spin-neutral 
electron pairs of momentum tt. The corresponding energy eigenstates are outside the 
Bethe-ansatz subspace. The remaining M^^ = [Na — 2Sc — 2Srj] ry-spinons are anti-bound 
within composite entities with total 77-spin zero (anti-bound state between a r7-spin-l/2 
?7-spinon with component 1/2 and another one with component —1/2). Again, there are 
as well Bethe-anstaz excited states containing = 1, 2, 3, ... ry-spin-neutral pairs of such 
anti-bound ry-spinons. Each pair involves two sites, doubly occupied and unoccupied by 
rotated electrons, respectively. 

The eigenvalue 2Sc of the operator that counts the number of rotated-electron 
singly occupied sites obeys the inequality 2Sc < N. As a simple example let us consider 
ground states with electronic density n = N/Na and spin density m = [N^ — N^l/Na in 
the ranges n G (0, 1) and m G (0,n), respectively. For such ground states one has that 
2Sc = N. Within the operator formulation of [Ij, those have Ms = [M^^ + M^^] = A^ 
spinous of which M^^ = 2Ss are unbound spinous and M^^ = 2Nsi are bound spinous 
inside A^^i = [A^/2 — Ss] spin-neutral two-spinon composite si fermions. Furthermore, 
such ground states have = 25^^ = [A^^ — ^] '^-spinous, Nc = N c fermions, and 

~ — N] c fermion holes. For them the number of electrons A^ equals that of 
rotated electrons that singly occupy sites and the number of rotated-electron doubly 
occupied sites vanishes. The ground-state Mg = N spinous refer to the A^ spin-1/2 
spins of the rotated electrons that singly occupy sites. The ground-state Nc = N c 
fermions describe the charge degrees of freedom of such rotated electrons. The ground- 
state = [Na — N] ry-spinons describe the ry-spin degrees of freedom of the [A^^ — 
sites unoccupied by rotated electrons. The ground-state A"^^ = [A^^ — A^] c fermion holes 
describe the c hidden U{1) symmetry degrees of freedom of the latter [A^^ — ^] sites 
ground-state occupancies. 

Consistent with the results briefly reported above, within the formulation of [Ij 
the 7y-spin degrees of freedom of the rotated-electron unoccupied sites with component 
1/2 and the rotated-electron doubly occupied sites with component —1/2 tend to be 
anti-bound, whereas the singly occupied sites by electrons of opposite spin projection 
tend to be bound. One expects therefore correlations that should decrease somewhat 
fast with distance between the members of each pair of ry-spinons or of spinous. 

The relevance of the correlations between doubly occupied sites and unoccupied 
sites has been suggested before by several authors. This can involve the introduction of 
an eflFective low-energy theory that contains a charge 2e bosonic mode [H EJ [9] , which 
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may be bound to a hole. The significance of short-range correlations between unoccupied 
and doubly occupied sites was, for instance, found in [TO]. The proposal of (anti-)bound 
states of doubly occupied sites and unoccupied sites in the repulsive Hubbard model at 
half-filling (Mott insulating phase) and of spins of opposite projections in the negative 
case (Luther-Emery phase), were recently justified by the existence of long range-order 
in a non-local order parameter pT]. 

The unitary transformation, V{U)^ is such that [i), H] = 0, where D = VDV~^ = 
VDV^ ^ with D the operator that counts the number of doubly occupied sites of the 
original electrons. At very large [/, the eigenstates of the Hamiltonian may be labelled 
by the eigenvalue of D and at finite U they may be labbelled by that of D. The rotated 
electrons (and in general any operator written in terms of the rotated electrons) can be 
obtained in the form c = VcV~^^ c^ = VcW~^ . The eigenstates of the Hamiltonian at 
any value of U may be generated as, 

\i^u) = V\iJu=c.) (1) 

which uniquely defines the unitary transformation [T]. Therefore any correlation 
function of two operators A and B satisfies, 

{i^u=oo\AB\iJu=oo) = (MVABV-'liJu) 

= {MAB\^u) (2) 

As a consequence any correlation function of operators involving rotated electrons at 
finite U may be obtained from the correlation function of the original electrons calculated 
at very large U . The form of the unitary transformation, V ^ is rather involved and 
contains infinite terms if expressed in terms of the original electron operators. However, 
to calculate their correlation functions it is enough to calculate the correlation function 
of the original electrons at large [/, which greatly simplifies the problem. (As a 
consequence, the correlation function of AB is independent of U .) It is the purpose 
of this work to calculate the real-space correlation functions of these objects. 

In this paper we use several methods to determine these correlation functions 
such as a mean-field theory based on the Zou- Anderson transformation [12j, a possible 
description in terms of an exact spin-charge-like separation of the original degrees of 
freedom introduced by Ostlund-Granath [13j and the density matrix renormalization 
group (DMRG) technique. 

The first issue considered in the following has to do with the definition of which 
correlations functions one wants to calculate. Those typically involve products of 
several electronic operators. Considering low energies, where a bosonization approach 
should apply, we expect that as the number of fields increases the absolute value of the 
correlation functions exponents should increase (considerably). In the metallic phase 
bosonization predicts a power-law decay with distance. (If the exponent is high then 
the extent of the correlation function should be very small). In the case of half-filling 
Umklapp scattering may change the behavior to an exponential decay. (However, if 
the exponent of the power-law decay is high the two behaviors will be to some extent 



Spinon and rj-spinon correlation functions 



5 



similar) . One of the aims of this work is to determine the exponents of these decays or 
their correlation lengths. 

2. Correlation functions 

The correlation function, 

Ci(r) = ((1 - nt(r)) (1 - n;(r)) n^{r = 0)n;(r = 0) (3) 
and its connected function, 

Ctir) = Ci(r) - ((1 - n^{r)) (1 - H{r))){n^{r = 0)n^{r = 0)) (4) 

contain information about the correlations between a unoccupied site (of the original 
electrons) at point r and a doubly occupied site at the origin. This is a charge correlation 
function. It is related (at large U) with the proposed (anti-)bound states of ry-spinons 
with opposite 77-spin projections. Those refer to the ry-spin degrees of freedom of pairs 
of rotated-electron doubly occupied and unoccupied sites. 
The correlation function, 

C2(r) = ((1 - nt(r)) n^(r)nt(r = 0) (1 - n^(r = 0))) (5) 

and its connected function, 

Clir) = C,ir) - ((1 - nt(r)) n^(r))(nt(r = 0) (1 - n^(r = 0))) (6) 

are related (at large U) with the proposed bound states of spinous with opposite spin 
projections. This is a spin-like correlation function. Even though the connection to the 
Bethe-ansatz states is through the rotated electrons, and only in the large U limit they 
are close to the original electrons, we will consider in this work the correlation functions 
at different values of U. 

We also calculate some mixed correlations where at point r we have for instance a 
doubly occupied site and at site r = a singly occupied site, or at site r a unoccupied 
site and at site r = a singly occupied site. That is we also calculate, for instance, 

C3(r) = ((1 - n^{r)) (1 - n^{r)) n^{r = 0) (1 - n^{r = 0))) (7) 

or 

C^(r) = {n^{r)n^{r)n^{r = 0) (1 - n^{r = 0))) (8) 

and the corresponding connected correlation functions. 

The operator that counts the number of doubly occupied sites may be written 
as Z^r- 4,t^^5t4,^^^,^' similarly for the unoccupied sites Z^^- <^r,t4,t^^5^4,^' singly 

occupied sites with spin up J2r 4,t^^,^^^r,^^^,t singly occupied sites with spin down 
J2r ^l,i^r,t^l^^^r,i' The cvaluatiou of the correlation functions depends on the band- filling. 
Away from half-filling (metallic phase) we expect from the bosonisation that both the 
charge and the spin correlation functions will decay with distance, as power laws. At 
half filling it is expected that the charge correlation functions become exponential-like. 
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due to the presence of a charge gap. In finite magnetic field the spin degrees of freedom 
wiU also develop a gap. 

The corresponding correlation functions are then expected to be of the form, 



for the charge correlation functions, if there is a charge gap such as at half-filling, where 
^ is the correlation length, and with a possible extra oscillating factor of the type (— 1)^. 
On the other hand, the spin correlation functions and the charge correlation functions 
in the metallic phase are expected to be of the form. 



also with a possible extra oscillating factor of the type (—1)^. 

In the half-filling case and in the limit of large U the spin part of the Hubbard 
model reduces to the spin-1/2 isotropic Heisenberg chain. Also the charge part is gapped. 
Previous studies for the charge correlations suggest that a ^ 1/2. The study of the spin- 
spin correlation functions at half-filling lead to some controversy about the presence of 
logarithmic corrections, but the presence of a logarithmic factor with exponent /5 = 0.5 
was confirmed and the above a exponent reading a ^ 1 [HI [151 [El [El [IE] . 

There are transformations proposed in the literature that lead to a similar 
decoupling of the electronic degrees of freedom. Examples are for instance given in 
[12] or in The main motivation was the study of either the large- [/ limit in the 

Hubbard or Anderson models [20], with the intent of controlling in a efficient way the 
projection to states where double occupancy is restricted (as in the t — J model), but 
considering a finite value of U instead of the extreme case of infinite [/, usually taken 
care of by a single slave boson [21j. Both representations introduce explicitly operators 
related to the four possible states associated with each site. Namely that a site may 
be unoccupied, singly ocuppied with a given spin projection or douby occupied. In 
the Kotliar and Ruckenstein procedure four bosonic operators are added, enlarging 
the operator space, that act as projectors on the original fermionic operators. In 
the Zou- Anderson transformation the original electron operators are replaced by two 
sets of two bosonic and fermionic operators that fulfill the projection. However, both 
representations lead to an enlargement of the physical Hilbert space and the extra 
unphysical states have to be projected out. We note, however, that the representation 
introduced by Zou and Anderson (ZA) has been used to explicitly obtain an exact 
solution of the Hubbard model in the large U limit in a much simpler way as compared 
to the Bethe ansatz [22]. Also, it has been used to study the stiffness of the one- 
dimensional Hubbard model in a way equivalent and alternative to the Bethe-ansatz 
solution [23j. 

2.1. Zou- Anderson transformation 



C{r) 



1 _r 

— e « 



(9) 



C{r) 




(10) 



The electron operators may be written as. 
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where the operators e^, d^, Si^a annihilate sites that are unoccupied, doubly-occupied and 
singly-occupied with an electron with spin a, respectively. Since the electron operators 
are fermionic we can either choose the operators e, as bosonic and the operators Si^a 
as fermionic, or vice-versa. In the original paper [12j the first choice was made (called 
slave-boson approach) and in [22] the second choice was taken (called slave-fermion 
approach). Most expressions are the same, formally, in either case. The difference 
arises when one integrates over degrees of freedom or in the mean-field approach when 
the Bose-Einstein or the Fermi-Dirac distributions appear. 

The enlargement of the degrees of freedom imposes the constraint, 

elei + dld, + Y,SlS^,, = l (12) 

cr 

at each site. (This is the completeness relation of the four possibilities: one site is 
either unoccupied, doubly-occupied, or is singly occupied by an electron with spin-up 
or down). The constraint is simply obtained imposing the anticommutation relation 
= 1 and considering either the slave-bosons or slave-fermions commutation 
or anticommutation relations. 

Let us consider the Hubbard model written as. 

The first is the hopping term between a site i and its neighbors distant by S (in general 
vectors in a rf-dimensional space) , U is the on-site repulsion and /x the chemical potential 
enforcing the band filling. 

In terms of the slave-bosons or slave-fermions the Hubbard Hamiltonian may be 
rewritten as, 

H = — t^^ ~ d^id]^^^ Sj^^Si^s,cr 




+ UY^ did, + /i ^ (eje, - did,) - fiNa (14) 

i i 

The ZA mapping reverses the role of the interacting and kinetic terms in the 
Hamiltonian. The interacting Hubbard term becomes quadratic in the ZA particles and 
the kinetic one is transformed into an interacting quartic term that couples particles 
along the lattice links. This is particularly useful to study the strongly interacting 
(large U) regime where the kinetic term is treated as a perturbation. The price of this 
transformation is the appearance of an on-site constraint, which assures exactly one 
particle per lattice site. In the mean field (MF) approach this translates to an on-site 
Lagrange multiplier. 

The problem to be solved involves the effective Hamiltonian, 

H = — t^^ (cicl^^ — didl^^) Sj^^Si^s,a 
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+ d\di + /X ^ {^A^i ~ A^i) — I^N 

i i 

+ E (^e. + did, + E S^S.^a - 1 j (15) 

where we have introduced at each site a Lagrange multipher, A^, to inforce the constraint. 
The transformation of the electron operators to the auxihary operators already embodies 
part of the classification of the Bethe-ansatz states (actually for the rotated electrons). 
It seems therefore natural to decouple the quartic terms in the Hamiltonian in such a 
way that the unoccupied and doubly-occupied sites are separated, on a first stage, from 
the singly-occupied sites. Also, we consider that the unoccupied sites and the doubly- 
occupied sites are paired on nearest-neighbor links. On the other hand, the spin states 
of the singly-occupied sites are paired into spin singlets. 

We then consider the mean-field Hamiltonian in the following form [24j . 

Hmf = ~* E { (^^ ~ Xs) Sl^Si^s.a 

+ (el+^e, - 4^gdi) xla - {xt - Xi) Xt} 



i,5 



+U J2 dUi + E (4ei - dUi) - nN 

i i 



(16) 



The quantities appearing in this Hamiltonian expression are defined as foUows, 

Xs = (el+s^i) 

xi ={dlA) 

= {eidi+s) 

^6 = {{Si^iSi^s,t - Si^-fSi^s^i)) (17) 

Besides considering hopping amplitudes we also introduce two pairing terms [24j, one 
between the unoccupied and doubly-occupied sites and another one between singly 
occupied sites with opposite spins. Note that one refers to a boson pairing and the 
other to a fermionic pairing. The choice of the mean-field parameters has in mind 
the possible bound states between the e and d operators and the and operators. 
Indeed, we intend to investigate the tendency to form these bound-states. The problem 
is now quadratic and may be diagonalized. The solution is briefiy reviewed in Appendix 
A. 
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Generically the phases found by solving the MF solutions for arbitrary band filling 
and energy are characterized as follows [24j: Phase (1) is conducting and characterized 
by X 7^ 0, A = 0. In it the spinous are gapless and the charge degrees of freedom exibit 
a gap of the order of the temperature, which closes at T = 0. This is the lowest free- 
energy phase. Within the mean-field approach it is that corresponding to the ground 
state. At finite energies other phases emerge [24], such as phase (2) x = 0, A 7^ 0, 
which is gapped for both degrees of freedom. Since it appears near x = it is tempting 
to identify it with an insulating antiferromagnet. Phase (3) x 7^ 0, A 7^ is a precursor 
of the superconductor. In it there exists spin-singlet formation but the charge motion 
is incoherent since no condensation is allowed. If one imposes ek=o = this phase 
splits into two sub-phases, analog to the pseudogap and superconducting phases in [25j. 
Phase (4) is an incoherent high-temperature phase where all correlations are zero. 

The hopping and pairing correlation functions between two sites at distance r from 
each other are given by, 

XF{r) = (4,1^0,1 + 4,-i^o,-i}q, 
XB{r) = (^d% - elco)^, 

AB{r) = {dreo + erdo)^. (18) 

They have been calculated before in Ref. [23] for the various phases. Although 
these correlation functions are not gauge invariant, they are useful to characterize 
the different phases. For phases (1) and (3) and both for the fermion and the boson 
hopping correlation functions, it was found close to half filling that the correlation length 
increases as the doping increases [24j . Particularly, the bosonic correlation function has 
a large correlation length. Analyzing the correlation length of A^, one clearly sees a 
long-range correlation in the high doping regime (x = 1 — n large), possibly precursor of 
Bose-condensation and superconductivity. In the low-doping region, both the bosonic 
and the fermionic correlation functions have a smaller range consistent with a spin 
gapped state. In this regime the two correlation functions have similar range, while at 
higher doping the charge correlation function has a much larger range compared to the 
spin correlation function. These correlation functions will be relevant if the system is in a 
deconfined phase. In a confined phase these correlation functions loose their significance, 
as the various degrees of freedom are confined within the real electrons. However, use 
of the Bethe-ansatz solution reveals that some fractionalization and rearrangement of 
the degrees of freedom occurs. 

We can now calculate the above mentioned correlation functions, Ci(r),C2(r), in 
the mean-field approach. 

Using the constraint we can write that, 

Ci{r) = ((1 - n^(r)) (1 - n;(r)) n^{r = 0)n^{r = 0)) = (ne(r)nrf(O)) 
C2{r) = ((1 - nt(r)) n^{r)n^{r = 0) (1 - n;(r = 0))) = (n,;(r)n,t(0)) 

(19) 
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Figure 1. Correlation functions (r), C2{r) for n = 0.5, U = 4, 10 and N = 100 in 
the phases 2 and 3. 
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Figure 2. Correlation functions (r), C2(r) for n = 0.7, U = 4, 10 and N = 100 in 
the phases 2 and 3. 



and the corresponding connected functions read, 

Ctir) = (ne(r)n,(0)) - (ne(r))(nd(0)) 

Clir) = (n4(r)n,t(0)) - (n4(r))Kt(0)) (20) 

where = e^e^ Ud = d)d^ = sjs^ and = s^^s^. Using their representations in 
momentum space in terms of the diagonalized operators we obtain the results derived 
in Appendix A. Those are shown in Figs. [T]and[2j 

At half filling we are in the insulating phase (2). We expect therefore a charge gap 
and an exponential decay of the correlation function Ci{r). This is indeed seen for the 
values of U = 4, 10, and independent ely of the system size, where the spatial extent 
refers basically to nearest-neighbors. The spin part is gapless, so that a larger range 
correlation function is expected, as shown in Fig. [3| As [/ increases, the magnitude of 
C2{r) grows for r = 1, but decreases faster with distance as compared to smaller values 
of U. Away from half filling we consider the densities n = 0.7,0.5. In the metallic 



Spinon and rj-spinon correlation functions 



11 



0.04 



0.03 - 



0.02- 



0.01- 



C^(r) 



n=l,Phase2 



— (ill i $ i (|i i $ i (|i — 




0.08 



0.06 



0.04 - 



0.02 



0.04 



0.03 



0.02 



0.01 



n=l, Phase 3 



2468 10 2468 10 

r r 




0.08 



0.06 



0.04 



0.02 



2468 10 2468 10 

r r 



Figure 3. Correlation functions (7f(r),C2(r) for n 
the phases 2 and 3. 



1,U = 4,10 and N = 100 in 



phase the charge correlation function has a much larger range, comparable to or larger 
than the spin counterpart. This qualitatively agrees with the results for the non-gauge 
invariant correlations. Far from half filling (quarter filling, n = 0.5) the spin correlations 
decrease fast with distance since we move far from the half-filled antiferromagnet. 

As stated above, there are two sorts of approximations within the present approach. 
The first is related to the enlargment of the physical Hilbert space and the necessity 
of introducing a constraint, to reduce the system to that space. The other sort of 
approximation has to do with the mean-field approach used. Moreover, the constraint 
is only implemented on average, as usual in slave-boson or slave-fermion approaches. 

The first difficulty has been overcome recently [13j , with the introduction of an exact 
transformation of the electron operators in terms of other operators that are related to 
the spin-charge separation of the model. The electron operators can be written as 
composites of charge-like and spin-like operators that do not give rise to any unphysical 
states, and thus avoids introducing any constraint. 



2.2. Ostlund-Granath transformation 

This transformation [13j introduces new operators called quasicharge and quasispin 
operators ql that obey, respectively, Fermi and Bose statistics, 

- (^_^^ ^, (21) 

(22) 



= 4 ,.(1 - n^^r) + {-lyc^^^n^,, 



— (c| ( 1)''C| ) c^^ 

= m 



These operators satisfy the algebra {c^,c^/} = S^y , {cj,c|!'} = , 



, 



i J 



iSrr' Y^k ^ijkQr- The clcctrou operators are expressed in terms of them as. 



(23) 
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The Ostlund-Granath representation involves other operators such as the quasicharge 
operator = clc^ and the local pseudospin operators pi = n^q^^ which are the 
generators of the SU (2) algebra that corresponds to "rotations" between the unoccupied 
and doubly occupied states [26] . 

The following results hold = 1 — 2n^g^, = (1 — ) ql and = [rir — 1)^ 
where s\ = | c^,r^a/5^/3,r5 with (7^ the Pauli matrices. The total z-component of 
pseudospin can therefore be seen to be half the number of doubly occupied sites minus 
the number of unoccupied sites, which is precisely the charge relative to half filling. The 
action of these operators onto the four-state basis is shown in Fig. [4} 

The Hubbard model can be rewritten in terms of the quasiparticle operators as, 

= t (To + Ti + T_i) + Uhu (24) 

with hu = iT^r ^l^ri 

^0 = ^ E (1 + 49r • + CC) (25) 

er,r' 

ri = ^E(-ir (1-49.- 9.0(44') 

and T_i = T/. 

Even though this transformation achieves some sort of exact spin-charge separation 
(actually quasispin and quasicharge), the Hamiltonian has a complicated structure. 
Although it has quart ic interacting terms, since they involve the quasispin operators, 
these are typically represented in terms of bilinear representations of fermionic or bosonic 
operators. Usually these representations enlarge the physical Hilbert space and one has 
to introduce constraints. One can use a Majorana fermion representation [27j but this 
leads to an Hamiltonian where the leading interacting term involves six operators and, 
therefore, the analytical treatment is rather complicated. One can also represent the 
quasispin operators using the Jordan- Wigner transformation. The transverse terms are 
linear in terms of new fermionic operators (the strings cancel out since only nearest- 
neighbor hoppings are considered) but the longitudinal term is again a bilinear in 
fermionic operators (which leads again to terms with six operators). An analytical 
treatment would need some approximation scheme, which is known to not yield good 
results in one-dimensional systems. 

It is interesting however to look at the behavior of the correlation functions of these 
operators. Moreover, it has been shown that the quasicharge operator is associated with 
a recently found hidden U{1) symmetry of the Hubbard model ^ (on any bipartite 
lattice). Together with the exact Bethe-ansatz solution of the ID problem, such a 
symmetry has lead to a deeper understanding of the physics of the model, including 
an understanding of the dressed scattering matrix structure [Ij. Specifically, we are 
interested in the correlation functions for the quasicharge and the quasispin operators 
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qz=-l/2 



Figure 4. Eigenvalues of the quasicharge, quasispin, pseudospin and spin for each of 
the four basis states referring to a doubly occupied site, singly occupied site with spin 
up and spin down and unoccupied site. The transitions generated by the off-diagonal 
operators are also shown. In the last panel we consider the transitions generated by 
Uf = I — Tie where f = c'^ ^ = c, similarly to the operators introduced in [28] (for the 
rotated electrons). 

of Ostlund and Granath written in terms of the original electron operators (and the 
corresponding connected correlation functions), 

C5(r) = (n,(r)n,(0)) = ((n(r) - if (n(0) - 1)') 

Ceir) = ((1 - ne(r)) (1 - n,iO))) = 1 - (n,(r)) - (nM) + K(r)ne(0)) 

C7(r) = {q.{r)qM) = ((^ " ^iW) " ^^(0))) 

(26) 




Figure 5. Connected correlation functions Ci (77-spinon), C2 (spinon), Cs and C4 
(spinon-77-spinon) . 

where n(r) = n^(r) + n^(r) and ria = cJ.Ccr. 

We may also rewrite the pseudospin correlation functions in terms of the original 
electron operators (and the corresponding connected correlation functions), 

Csir) = {p.{r)pM) = (Hr) - if - n,(r)) (n(0) - if - n,(0))) 

(27) 

The spin correlation function can be evaluated in the usual way, 

C9(r) = {s,{r)s,{0)) (28) 

A direct solution of these correlation functions in terms of the Ostlund and Granath 
Hamiltonian is complicated, since even the mean-field approach is complex. Therefore, 
we have used a DMRG method to calculate these and other correlation functions. 

2.3. DMRG calculations: Correlation functions 

For simplicity, here we limit ourselves to the half-filling case. Using the DMRG method 
(briefiy reviewed in Appendix B) we have calculated the correlation functions indicated 
above. We consider correlation functions for the original electron operators as a function 
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of U. We emphasize that in the hmit of very large U these equal the correlation functions 
of the rotated electrons for [/ > 0. 

The charge correlation functions typically decay fast as the distance r between the 
two operators grows. On the other hand, the spin related correlation functions oscillate 
by a factor of the type (—1)^ and decay slowly with distance with a power law behavior, 
like in the spin- 1/2 isotropic Heisenberg model. 

In Fig. ([5]) we present results for the correlations functions Ci, C2, C3 and C4. 
The correlation function Ci{r) decays fast with distance, which shows that the doubly 
occupied sites and the unoccupied sites are tightly correlated in the half-filled phase. As 
U grows the spatial extent is strongly reduced and in the large U limit the correlation 
is basically extended to the nearest-neighbors. Since in the very large U limit a doubly 
occupied site costs an infinite energy, the correlation function basically vanishes. The 
very large U limit corresponds to the correlation function of the rotated electrons for the 
whole U > range. Hence this implies a very short range and a vanishing correlation 
function in the strictly infinite limit. 

The correlation functions for the doubly occupied sites and a singly occupied site 
or a unoccupied site and a singly occupied site are negative. This is indicative of an 
anti-correlation, as expected. These correlations also decay very fast with distance. 
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On the other hand, the correlation function C2{r) for a singly-occupied site with 
spin up and a singly occupied site with opposite spin projection, is somewhat similar 
to the longitudinal spin correlation function. It oscillates with distance and decays 
slowly with r, which indicates a long-range correlation. The influence of the Hubbard 
interaction is smaller than for the charge correlation functions. For large U the short 
range values of the correlation function increases, consistently with a more pronounced 
spin character of the excitations of the half-flUed Hubbard model at large interactions. 
However, as U becomes very large (for instance comparing U = 100 with U = 20) the 
decay with distance is faster. 

In Fig. ([g]) we show results for the correlation functions for the quasicharge, 
quasispin and pseudospin operators introduced in Ref. [13j. The charge correlation 
functions decay fast in a way similar to Ci(r), speciflcally the quasicharge and the 
pseudospin correlation functions. The quasicharge shows a correlation while the 
pseudospin shows an anti-correlation. The two spin correlation functions have a slowly 
decaying oscillating behavior. Analysis of both the quasispin and the spin correlation 
functions reveals that the nearest-neighbor is anti-correlated. On the other hand, 
the spinon correlation function C2{r) behavior shows that the nearest-neighbors are 
positively correlated, since it is an occupation number correlation function. The other 
correlation functions respect the opposite spin projections. 

The charge correlation functions can be fltted with an expression of the form given 
in Eq. ([9]). Both the decay length ^ and the exponent a are indicative of a stronger 
decay as the coupling grows. The results are shown in the Table. 

[/ = 4 U = 10 U = 20 U = 100 
a 2.16 4.84 6.43 11.17 
^ 2.03 1.39 1. 1. 

The spin correlation functions can be fltted to an expression of the form provided in 



Eq. (10). The results for the various spin related correlation functions show that their 
decay is very similar to the decay of the spin correlation function Cg{r). In the inflnite 
U limit these correlation functions tend to the corresponding correlation functions of 
the spin- 1/2 isotropic Heisenberg model. 



2.4' Eigenstates of the reduced density matrix 

The correlation functions we have calculated give information on the correlations 
between the various entities discussed above. A more direct approach is obtained 
studying the eigenstates and eigenvalues of the reduced density matrix of two sites 
in the chain. Considering a chain of sites we may single out two sites distant by r 
lattice units. The full density matrix of the chain can then be written as, 

Piv = |V')(V'| (29) 

where is the ground state of the system. The ground state may be written as the 
direct product of the states at each site. These can be written in terms of a basis with 
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four states, namely |0) = |0), | t)? I 4-) J t-Dj referring to the four possibilities that each 
site is either unoccupied, occupied by a particle of spin up, a particle of spin down or 
doubly occupied, respectively. The full density matrix is a 4^ x 4^ matrix. The ground 
state may be obtained for instance considering exact diagonalization of small systems. 
We have used Lanczos method to obtain the ground state expressed in this basis. We 
considered a system of size = 14. 

Information about the correlation between two points on the lattice may be obtained 
considering a reduced density matrix by integrating N — 2 sites. One of the sites may 
be located at point r = on the lattice and the other may be located at site r. The 
reduced density matrix is then obtained as, 

P2 = TrA^-2|^)(^| (30) 

This is a 16 X 16 matrix that can be diagonalized for different values of U and different 
values of r. The eigenvalues give the probabilities to find the two sites in a given 
correlated state characterized by the corresponding eigenstate. We consider as before 
U = 4,10,20,100 and r = l,---,7. In addition, we consider half-filling and zero 
magnetization, which implies seven electrons with spin up and seven electrons with 
spin down in a chain with fourteen sites. 



eigenstate |0r=o;0r) = |^o,t^r,t^o,;^a) 



A a (|1001) + |0110)) + /3 (|1010) + |0101)) 

B ^(10110) -11001)) 

C |1100) 

D joOll) 

E ^(11011) + 10111)) 

F _!= (11110) + 11101)) 

G ^(|1000) + 10100)) 

H ^ (lOOlO) + 10001)) 

/ ^(11011) -10111)) 

J _!= (11110) -11101)) 

K ^ (|0010) - |0001)) 

M |1111) 

TV 1 0000) 

O _!= (11010) -10101)) 

P 7 (|1001) + |0110)) - 5 (|1010) + |0101)) 



The structure of the normalized eigenstates is illustrated in the table. Those are 
linear combinations of a few states of the basis for sites i, j represented as {ni^^nj^^ni^^rij^i) 
where rii = 0, 1 and means unoccupied and 1 occupied. Since we are considering 
half filling, there is particle-hole symmetry. Moreover, zero magnetization implies a 
symmetry between up and down spins. The eigenstate with the highest eigenvalue 
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Figure 7. Eigenvalues of the reduced density matrix, p2, for ^7 = 10 (left panel) and 
U = 20 (right panel) as a function of distance r, for some eigenstates. 



is the state A represented in the table. This means it has the highest weight in the 
correlations between the two sites. 

State A has four components with relative weights a, /5. The coefficient a measures 
the contribution due to spins up and down in the two sites and the coefficient /3 gives 
the contribution of a double occupied site and a unoccupied site at r = and r. In 
Fig. [7| we compare the weights of some of the eigenstates as a function of the distance 
r for [/ = 10, 20. This state is spatially symmetric (exchanging the two sites) and anti- 
symmetric in spin space. It corresponds to a spin singlet. The part of the state with 
weight a is a spin singlet between sites r = and r and the state with weight /3 is a 
local spin singlet at either site r = or site r and a 77-spin triplet by exchanging the 
doubly occupied site with the unoccupied site. At half filling we expect that for large U 
the state A has a large weight, since it is associated with the formation of spin singlets. 

State O is associated with an antisymmetric pairing of a doubly occupied site with 
a unoccupied site. It has a local spin singlet at the doubly occupied site. This state has 
very small weight whereas state P is the lowest-weight state. It is the counterpart of 
state A. 

States B and C (and its degenerate state D) also have a large weight. These states 
are spin triplets and antisymmetric in the space coordinates. Interestingly states where 
the spins at points r = and r are parallel have equal or larger weight. For odd sites, 
states B^C^D are degenerate but for even sites it is more favorable for the spins to be 
parallel. This is indicative of the long-range antiferromagnetic correlations in the large 
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0-71 1 




r 



Figure 8. Relative weight of the spin up-spin - down-spin subspace a and doubly 
occupied-site - unoccupied-site subspace /3 in the eigenstate with highest eigenvalue, 
A, for the matrix p2. 

U limit at half- filling. 

States E and / are representative of mixtures between states involving singly 
occupied sites and doubly or unoccupied sites, either in antisymmetric or symmetric 
combinations, and have smaller weights. 

Entanglement between doubly occupied sites and unoccupied sites enters indirectly 
via state ^4, since state O has very low weight. The lowest probability state P also has 
a mixture of the same four states, as well with a very small weight. 

Even though the doubly occupied-site - unoccupied-site entanglement contributes 
to state the relative weight decreases fast with distance and interaction strength. This 
is illustrated in Fig. |8| In the left panel we show the relative weight of the spin up-spin 

- down-spin subspace a and in the right panel we show that of the doubly-occupied-site 

- unoccupied-site subspace /3. Note that the dominant contribution comes from the spin 
up-spin - down-spin subspace a > 13. (As a side remark, the eigenstate with the lowest 
eigenvalue, P, has 7 < S.) Comparing the various values of the interaction strength [/, 
we see that /3 decreases fast. The same happens as r increases, which is consistent with 
the short range of the correlations between a doubly occupied site and a unoccupied 
site. For large U the weight is mostly contained in the spin part, in a spin-singlet state. 
These results are consistent and clarify the previous results that the charge correlations 
are very short range, particularly as U grows. 



Spinon and rj-spinon correlation functions 



20 



Tracing out all states except those with unoccupied or doubly occupied sites, we find 
a 4 X 4 matrix that stores direct information on the correlations between a unoccupied 
site and a doubly occupied site. This density matrix may be defined as, 

Ps = Tr'p2 (31) 

where the trace is over the singly occupied sites. 

The eigenstates of this reduced density matrix are of the form. 



eigenstate |(/)r=o;0r) = |^o,t^r,t^o,;^a) 

/ ^(11010) + 10101)) 

// |1111) 

/// 1 0000) 

IV -1= (11010) -10101)) 

The results for the eigenvalues for U = 10, 20 of this reduced density matrix are shown 
in Fig. |9j The eigenvalue coresponding to state / that mixes a doubly occupied site 
with a unoccupied site is the largest at small distances and as distance increases all 
states become equally probable. The decrease of the relative weight of state / increases 
with [/, consistent with previous results. This state is a ry-spin triplet (note that this 
is consistent with the structure in Eq. 16). The lowest-weight state is a 77-spin singlet. 
The degenerate states // and /// correspond to two doubly occupied sites and two 
unoccupied sites. 

For the subspace spanned by {|0 >, | i>} x{|0 >, | i>} defined on two sites, the 
entanglement can be measured by the concept of concurrence [29j . While the concurence 
is still a combination of some correlation functions, it is diflFerent from the traditional 
density-density or other types of correlation function. The entanglement results from 
the linear superposition principle of the quantum mechanics and is absent in the classical 
physics. Therefore, it is usually regarded as a kind of pure quantum correlation. 

The reduced density matrix can be written as. 



P3 



/ u 








\ 





Wi 


z 








z* 







I 








V ) 



(32) 



The concurrence, as the measure of the entanglement, can be calculated as, 

C = 2max[0, \z\ - ^] (33) 

The results for the concurrence show that both for U = 10 and U = 20 the correlations 
only extend to nearest-neighbors. The concurrence for higher values of the distance 
between the two sites r > 1 vanishes. For nearest neighbors the concurrence takes the 
values C = 0.9577 and C = 0.9884 for U = 10 and U = 20, respectively, showing the 
large increase for large U and that the correlation is strong, since the concurrence is 
close to 1. 
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U=10 U=20 




Figure 9. Eigenvalues of the reduced density matrix pa for = 10 (left panel) and 
U = 20 (right panel), as a function of distance r. 

3. Summary 

In this paper we have studied the correlation functions of basic entities of the one- 
dimensional Hubbard model using various methods. Previous analysis of the exact 
solution via the Bethe ansatz suggests the importance of correlations between doubly 
occupied and unoccupied sites and sites singly occupied with spin up and spin down 
electrons. These correlations have also been suggested by other treatments, as mentioned 
in the text. 

The relevance of these operators is also stressed by their connection in terms of 
so-called rotated electrons with basic entities of the exact solution such as spinous 
and 7y-spinons. We have, therefore, calculated various correlation functions using an 
approximate mean-field solution within the Zou- Anderson transformation and with the 
introduction of non-local bond variables, in a way that is a reminder of the spinon and 
77-spinon bound and anti-bound states, respectively. Furthermore, we have used exact 
DMRG calculations of the same quantities and supplemented those with calculations of 
correlation functions for the operators introduced by Ostlund and Granath, and further 
developed by one of the authors, which allow an exact charge-spin separation in a way 
that reminds the exact separation in low-dimensional systems. 

We concluded that the charge-like correlation functions are typically very short 
ranged in the case of half filling, due to the charge gap. As the interaction strength 
increases, the correlations become virtually nearest-neighbor like. The spin-like 
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correlation functions are however more extended in a way similar to the spin-spin 
correlation functions of the Hubbard model or, in the large U limit, those of the spin- 
1/2 isotropic Heisenberg model. Even though the correlation functions calculated here 
are different from a standard longitudinal spin correlation function, their decay with 
distance is similar. 

Further insight onto the correlations between two sites in a chain was obtained 
calculating the eigenvalues and eigenstates of the two-site reduced density matrix using 
exact diagonalization of a small system, also at half filling. The role of the spin and 
charge contributions to the entanglement between the two sites was clarified. As shown 
by the other methods, the correlations between doubly occupied sites and unoccupied 
sites are very short ranged, as evidenced by the concurrence which extends only to 
nearest neighbors. The eigenstates of the reduced density matrix corresponding to 
singly occupied sites with spin up or spin down have longer range. 

The results using the mean-field approach were also extended to cases in the metallic 
phase (away from half-filling). The absence of the charge gap leads to correlations 
between a doubly occupied site and a unoccupied site that have larger range. The 
mean-field treatment also allows the study of higher energy phases. It was shown that 
for some of these phases the charge correlations have comparable ranges to those of the 
spin correlations. 
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Appendix A. Mean-field solution using the Zou, Anderson transformation 

Here we briefly review the mean-fleld solution of the Hubbard model in terms of the 
link variables introduced in Eq. (16). The bosonic and fermionic parts decouple. The 
mean-fleld Hamiltonian may be written as, 

HMF = H,d + Hss + C (A.l) 

where 

1,5 

+ d\di {U - ^l + Xi)^Y. (a* + \) (A.2) 
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+ >^^SlS,^. (A.3) 

C-t j:{{xl-xt) Xl + ^sA*s + m^s} 

-fj,N-J2\ (A.4) 

i 

Defining tiie Fourier transforms of the operators as, 

er = ^E^"''-ek (A.5) 



and 



xt = E^-'''xf (A.6) 
s 

where A — e,d, S, and similarly for $fc and A^, we can write in momentum space that, 
^ed = E I (-^ExL + + Aoj elck +(^t^xla + U-^i + Aoj 44 
+ {-tAlekd-k-tAkdlelf,)} (A.7) 



That is, 



where 



Hed = E {^fee^efc + Dkdldk - tA*efcd_fc - tAkdlel^} (A.8) 
^fe = -^Exf,a + /« + Ao 

(J 

= ^ExL + t^-/" + Ao 



and 



Hss = E |E (-^(x^ - X^) + Ao) 

- t (^it^l,,^ - 4,4^1.,t) - m {Sk,iS.,,, - S,,,S.,,^)} (A.9) 
One tlien arrives to, 

k,a 

- tmSk,iS.k,t - Sk,^S.k,i)} (A.IO) 

wliere 

h = -t{xl-xt) + ^o (A.ll) 
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We have assumed that the Lagrange multipher is uniform and only the k — component 
is non-vanishing. For later purposes we define, 

s 

We are left with the diagonalization of two quadratic Hamiltonians. This can be 
done in a standard way performing Bogoliubov-Valatin transformations. Defining, 

d_k = Ukl3-k + vial (A.12) 
and similarly for the spin part, 

Si-k = Uk^-k - vl&l (A. 13) 

where u\ — vl — l for bosons and u\ + vl — l for fermions, and eliminating off-diagonal 
terms in the quasiparticle operators, we obtain that the diagonalized Hamiltonians take 
the forms. 



where 



Hss = E {^k (4«fc + ^-J-k) + [ek - h]} (A.14) 

k 

et = ek + 1/2{E, - D^) 

el = ek- ll2{Ek - D^) (A.15) 



with 



2 



Note that Ek + Dk = U + 2Xo. Here, 

h = \/^+4t^ (A.17) 



Also, we find that, 

UkVk 



Uk 



2ek 

1 f Ek + Dk \ 

2 1 2ek J 



vl (A.18) 



2 I 2efc 



and 



UkVk = — 



ul — 
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The solution of the problem involves the self-consistent calculation of the averages 
appearing in these equations. Any average can now be calculated using the Bogoliubov- 
Valatin transformations. 

The mean-field equations are given by, 



= Jf^'yk-k' {uk'Vk' (1 - 2/(efc/))} 



(A.20) 



« k' 



1 



= TT l]7fc-fc' { + Uk'Vk' (1 + 5k',k„Nsna + (1 - Sk',kJfB{4')) 
-'^s k' 

+ Uk'Vk' (Sk',kpNsni3 + {1- Sk',k^)fB{4'))} (A-21) 
Xla = ^ {4/(40 + vl (1 - fiik'))} (A.22) 



k' 



Xk 



^J2^k-k'{ + ul,{Sk',kaNsnc, + {'i--Sk',kJfB{et')) 

^ k' 

+ vl (l + 5k',k,Nsnp + (1 - 4',fc,)/B(4)) } (A.23) 

Xk = jrY.^k-k' { + ul {5k',k^Nsnp + (1 - 5k',kp)fB{ek')) 



1 



+ vl (1 + 5k',kc,Nsna + (1 - Sk',kJfB{et> 



(A.24) 



1 = ]^ E { + «fc (4,fc„A^«n, + (1 - Sk,kJfB{et)) 

^ k 

+ vl[l + 5k,k,N,np + (1 - 5k,k,)fB{el)) 

+ (Sk^^^Nsfip + (1 - 4,/c,)/B(ef )) 
+ i;,^ (1 + 5k,k^N,n^ + (1 - 4,.J/B(e^)) 
+ 2ulf{ek)-r2vl{l-f{ek))} 



n 



(A.25) 
(A.26) 



k-^k(y 



k^ki3 



where / is the Fermi-Dirac distribution and the Bose-Einstein distribution. 

The correlation functions may be calculated at the mean-field level. We obtain 

that, 

^iW = ^EE 

^^a ki k2 

\uk,\'\ukf + e-'^'^-'^^^ulvk,UkA^] fBK)fB{ei) 
\ukA'\vk,\' + e-'^'^-'^Klvk,UkA,] fB{el){l + Sb^)) 
\vk, ?\uk^' + e-^^^^-^->ulvk,UkA} + fBi^i)) 



+ 

+ 
+ 



\Vkr?\Vk, 



? + e''^'^-'^>ul^Vk,UkA^ 



+ + (A.27) 
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and 



Similarly, 



and 



k 



C2{r) 



/V2 



+ 
+ 
+ 



J(ki-k2)r~* 



(i-fKm-fieO) 



1 



T.\\^k\'M) + \v,\'{l-f{et)) 



Mr)) = ^ E [\^k\'f{et) + \vk\\l - f{el)) 

k 



(A.28) 



(A.29) 



(A.30) 



Appendix B. DMRG method 

In this Appendix a brief description of the DMRG method used in our studies is 
presented. The density matrix renormalization group (DMRG) [30] algorithm is an 
accurate method in dealing with quasi-one-dimensional system. It provides a criterion 
to find which states to keep and which to discard. Therefore, it can deal with relative 
large-size system with high accuracy. In the following, we give a brief introduction about 
this method and its application to our calculations. 

In this method one constructs a superblock composed of the original system block 
and the environment block, usually the refiection of the system block. The reduced 
density matrix for the system block is defined as, 

3 

where ^0 is a state of the superblock. Usually it is chosen as the ground state 
corresponding to the zero temperature. Moreover, here \i) and \j) label the states 
of the system and the environment blocks, respectively. For any system block operator 
A, we have, 



{A) = Y^^n'Pn' = Tr (pA) = T^^uj^ {u^ \A\ u^) . 



(B.2) 



where oja and are the eigenvalues and eigenstates of the reduced matrix p. The 
significance of the state u"^ can then be determined by uja- For a certain a, if c^c, is very 
small, its contribution to (^4) is also very small. Then its corresponding state can be 
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discarded. Using this method, some states can be discarded during the growth of the 
system size and the size of the Hamihtonian of the system to be calculated is therefore 
reduced. 

There are two basic DMRG algorithms — the infinite system and the finite system 
algorithms. For the infinite system case, the main process is as follows. We first choose 
a small-size system that can be exactly diagonalized, e.g. L = 4 sites is taken in our 
calculations, as the superblock. Then use its ground state to form the reduced density 
matrix p of the system block. The m highest eigenvectors of p are kept to renormalize 
the Hamihtonian of the system block S and the corresponding operators. We add then 
two new sites and use these renormalized Hamiltonian and operators of S to form a 
new superblock. By repeating these steps, the system size grows but the size of the 
Hamiltonian of the superblock keeps on a suitable size. 

The finite-system case is based on the infinite system case by sweeping the 
superblock to reach a higher accuracy. For further detailes on the process, see Ref. [30j. 

In our DMRG calculations we have used the finite-system DMRG algorithm. Three 
sweeps have been taken to increase the accuracy. The numerical calculations were 
performed on finite chains, up to 100 lattice sites, using the open boundary condition. 
Two sites were added in each step. For accuracy, the largest kept-state number reached 
130 and the truncation error is less than 10"''. To avoid the infiuence of the edge effect, 
the r = point was chosen in the middle of the chain. 
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